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I. INTRODUCTION 

The knowledge of electric dipole polarizability, a, of 
atoms and ions is required in different areas of physics 
and chemistry. It is the lowest order linear response prop- 
erty and relevant to a wide range of physical phenomena 
ranging from the microscopic to the macroscopic proper- 
ties. To mention a few macroscopic properties, the dielec- 
tric constant and refractive index of gas are among the 
important ones. In the case of microscopic properties, 
the parity non-conservation in atoms p], optical atomic 
clocks m [3] and physics with the condensates of dilute 
atomic gases [1H5] are of current interest. For accurate 
theoretical calculation of a, a precise treatment of the 
electron correlation effects is very important. In the past, 
a wide range of atomic many body theories were used to 
calculate a. The recent review by mitroy et al [7] gives a 
detailed overview of the atomic and ionic polarizabilities. 

In the present work we use the PRCC theory to cal- 
culate the a of the noble gas atoms. It is a theory we 
have developed to incorporate a perturbation in the con- 
ventional relativistic coupled-cluster (RCC) theory. In 
general, the coupled-cluster theory (CCT) [8l[9] is one of 
the most elegant many body theory which takes into ac- 
count the electron correlation to all order. The details of 
the CCT and different variants are described in a recent 
review [10]. The theory has been widely used for atomic 
[TTHTi] . molecular [TS], nuclear [Hj and condensed matter 
physics [T7] calculations. The PRCC theory is different 
from the previous RCC based theories in a number of 
ways. The most important one is the representation of 
the cluster operators in the PRCC theory is very differ- 
ent and it has the advantage of incorporating multiple 
perturbations of different ranks in the electronic sector. 
One basic technical advantage of PRCC is, it does away 
with the summation over intermediate states in the first 
order time-independent perturbation theory. In stead, 
the summation over all the possible intermediate states 
within the chosen basis set is subsumed in the perturbed 
cluster amplitudes. 

For the calculations we use the no virtual pair Dirac- 



Coulomb-Breit Hamiltonian. However, to assess the im- 
portance of Breit interaction we also carry out another 
series of calculations with the no virtual pair Dirac- 
Coulomb Hamiltonian. We isolate the changes aris- 
ing from the Breit interaction by comparing the results 
from the two sets of the calculations. In contrast, till 
date, majority of the theoretical calculations of a have 
been done with the Dirac-Coulomb Hamiltonian. We 
have chosen the noble gas atoms to study as these sys- 
tems are ideal to test the closed-shell PRCC theory. In 
previous works, a of the noble gas atoms were calcu- 
lated in the framework of many body perturbation the- 
ory HH], non-relativistic CCSDT [11] and RCC single, 
double and triple (RCCSDT) excitation approximation 
[5U] . In the later work, using RCCSDT, the third or- 
der Douglas-Kroll method [21] was used. It is an al- 
ternative of the Foldy-Wouthuysen (FW) transformation 
and a quasi- relativistic treatment. For the single particle 
wave functions, we use the kinetically balanced Gaussian 
type Dirac-Hartree-Fock orbitals. The results from our 
PRCC theory calculations are in good agreement with 
the experimental data and consistent with previous cal- 
culations. 

The paper is organized as follows. In the Section. 
II, for completeness and easy reference we briefly de- 
scribe the RCC theory with the Breit interaction. In 
section. HI we introduce the PRCC theory and pro- 
vide a detailed description of the tensor structure of the 
PRCC operators. In section. HI B we give the an- 
alytical structure of the PRCC equations. In section. 
Ill C we present a diagrammatic and algebraic descrip- 
tion of the the non-linear terms in the PRCC theory. 
In section. IV we introduce the formal expression of 
the dipole polarizability and its representation in the 
PRCC theory. In the subsequent sections we describe 
the calculational part, and present the results and dis- 
cussions. We then end with conclusions. All the results 
presented in this work and related calculations are in 
atomic units ( h = rUe = e = 47reo = 1). In this sys- 
tem of units the velocity of Hght is , the inverse of 
fine structure constant. For which we use the value of 
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137.035 999 074 ^H]. 



II. RCC THEORY 

For the high-Z atoms and ions, the Dirac-Coulonib- 
Breit HamUtonian, denoted by H^'~^^, is appropriate to 
include the relativistic effects. However, there are com- 
phcations associated with the negative energy continuum 
states of H^'^^. These lead to variational collapse and 
continuum dissolution [23j . A formal approach to avoid 
these complications is to use the no-virtual-pair approx- 
imation. In this approximation, for a neutral atom of N 
electrons [24 
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where a and /3 are the Dirac matrices, A+-(_ is an oper- 
ator which projects to the positive energy solutions and 
VAr(ri) is the nuclear potential. Sandwiching the Hamil- 
tonian with A++ ensures that the effects of the negative 
energy continuum states are neglected in the calculations. 
The last two terms, l/j'ij and g^{rij) are the Coulomb 
and Breit interactions, respectively. The later, Breit in- 
teraction, represents the inter-electron magnetic interac- 
tions and is given by 
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The Hamiltonian satisfies the eigen-value equation 
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where, is the exact atomic state and Ei is the energy 
of the atomic state. In CCT the exact atomic state is 
given by the ansatz 
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where is the reference state wave-function and r*-*^-* 
is the unperturbed cluster operator. In case of closed- 
shell atom the model space of the ground state consist 
of a single Slater determinant, |$o) ■ For an A^-electron 

closed-shell atom T^^^ — Y^^=i where i is the order 
of excitation. In the coupled-cluster single and double 
(CCSD) excitation approximation. 
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The CCSD is a good starting point for structure and 
properties calculations of closed-shell atoms and ions. In 
the second quantized representation 
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where t^^^ are cluster amplitudes, a\ {ai) are single parti- 
cle creation (annihilation) operators and abc. . . {pgr . . .) 
represent core (virtual) single particle states or orbitals. 
The eigenvalue equation of the closed-shell ground state 
in CCT is 



(7) 



Following similar procedure, the CC eigenvalue equation 
of closed-shell excited states may be defined as well. 



III. PRCC THEORY 

To incorporate an additional interaction Hamiltonian 
^^int perturbatively, we introduce the perturbed coupled- 
cluster operator T*^-*-^. This means, 7?int is applied once 
and residual Coulomb interaction to all order in all pos- 
sible sequences. In general, T^-*-) is a tensor operator 
and the multipole structure depends on the properties of 
i?int. The properties and values of T^-"-' indicate the ef- 
fect of Hint to the atomic state. With the perturbation, 
the modified eigenvalue equation is 
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where A is the perturbation parameter. Consider the case 
where Hint represents the interaction with an external 
static electric field E. The interaction Hamiltonian is 
then Hint = — J2i i"i • E = D • E, where D is the many 
electron electric dipole operator. The perturbed atomic 
state in PRCC theory is 
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This approach has the advantage of taking into account 
the effect of multiple perturbations systematically. Other 
than E, Hint could be one of the interactions internal 
to the atom like Breit interaction, hyperfine interaction, 
etc. For the present work, we examine T^^^ arising from 
E which is parity odd and vector in the electronic space. 





FIG. 1. Diagrammatic representation of T-,^ and T!, 



A. Tensor structure of PRCC operator 

For the present case, E as the perturbation, we can 
write the perturbed single excitation cluster operator as 



(1) _ 



TPCi{r)alaa. 



(10) 
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Note that T^"'^'' is a vector operator in the electronic space 
and the C-tensor Ci(f) represents the vector nature op- 
erator. The key difference of T^"'^-' from T^^^ is la + Ip 
must be odd, in other words (— l)'»+'p = —1. Here, la 
(Ip) is the orbital angular momentum of the core (vir- 
tual) orbital a (p). Diagrammatically, the T^^-* operator 
is represented as shown in Fig. [ija). It is similar to the 

conventional representation of T^^"^ but the interaction 
line is replaced by a wavy line. 

The tensor structure of T2^-' , on the other hand, has 
additional complications as it consists of two vertices. 
After due consideration of the iJi^t and T*^'''' multipole 
structure, it is represented as 
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Like in T^^'', are the C-tensor operators. Here, two 
C-tensor operators of rank I and k are coupled to a rank 
one tensor operator, Tj^-*. At the two vertices, the orbital 
angular momenta must satisfy the triangular conditions 
|ja - jp\ ^ I ^ (ja + 3p) and |jb - jq\ ^ k ^ {jb + jg). 
In addition, the two tensor operators must be such that 
1^ — fc| ^ 1 ^ (Z -I- k). These selection rules arise from the 
triangular conditions at the vertices. The other selection 
rule follows from the parity condition. For the present 
case as Hi^t is parity odd (— — _ (_i)('6+'<!). 

Only then the matrix element {pq\T'^^\ab) is nonzero. 
More details on the first principle analysis of the tensor 
structure, based on many-body perturbation theory, is 
given in Ref. 1251 . The diagrammatic representation of 

Tj^"* is as shown in Fig. [lib), where the vertical line 
on the interaction line is to represent the rank of the 
operator. Further more, this representation, at a later 
stage, simplifies the angular integration using diagrams. 



B. PRCC equations 

The eigenvalue equation with the perturbed Hamilto- 
nian is 

(12) 

When Hint is parity odd, like in the present case, there 
is no first order perturbative correction to energy Eq = 
Eq. In the CCSD approximation we define the perturbed 
cluster operator T'^^ as 
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Using this, the PRCC equations are derived from Eq. 
(12). The derivation involves several operator contrac- 



tions and these are more transparent with the normal 
ordered Hamiltonian iJ^^B ^ ^dcb _ ($^|iyDCB|$.^ 

The eigenvalue equation then assumes the form 
[Hn 



,^XH.J l*n) = \En - (<i>n|i/°^^l<i>n)] l*n 



.(14) 

A more convenient form of the eigenvalue equation is 

+ AiJint) |*o) - AEol^o), (15) 



where, AEq = Eq - {^o\H^'^^\^a) is the ground state 
correlation energy. Following the definition in Eq. (13 1, 
the PRCC eigen-value equation is 
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Applying e 



from the left, we get 



where ijDCB ^ 
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is the similarity 



transformed Hamiltonian. Using the Campbell-Baker- 
Hausdorff expansion 
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The commutations represent contractions and as iJ^CB 
consist of one- and two-body interactions, the expansion 



terminates at the fourth order. Multiply Eq. (17) from 



left by e ^'^ and consider terms linear in A, we get the 
PRCC equation 
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Here, the similarity transformed interaction Hamiltonian 
Hint terminates at second order as Hint is a one-body 



interaction Hamiltonian. 

Expanding H^'^^ and -ffint, the PRCC equation as- 
sumes the form 
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Here after, for simplicity, we drop E from the equations 
and for compact notation, we use H^ to denote 
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The cluster equations of T^"'^'' are obtained after project- 
ing the equation on singly excited states These 
excitation states, however, must be opposite in parity to 



|<I>o)- The T2"'^-' equations are obtain in a similar way af- 
ter projecting on the doubly excited states (^^bl- After 
the application of Wick's theorem, the cluster equations 
are 
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Where AB represents all possible contractions between 
the two operators A and B. The Eq. (21 ) and (22 1 form 



a set of coupled nonlinear algebraic equations. However, 
T'^^^ are solved first as these are independent of T^^^ the 
PRCC equations are then reduced to coupled linear al- 
gebraic equations. An approximation which incorporates 
all the important many-body effects like random phase 
approximation (RPA) is the linearized PRCC (LPRCC). 
In this approximation, only the terms linear in T, equiv- 
alent to retaining only H^^'^'> and r^^°'> in Eq. (21) 
and (22), are considered in the equations. Here after we 



use T as the general representation of the both T'^'^' and 
T^^' operators. 



side of Eq. (21) and (22) expands to 
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Here, T>t[°^ and BTl^'T^"' are nonzero only for T 

and T2^\ respectively. Each of the terms, after contrac- 
tion, generate several topologically unique Goldstone dia- 
grams. At this stage, the diagrammatic treatment is the 
preferred mode of further analysis and calculation. It 
simplifies the calculations and is well suited to represent 
contractions between the operators. 



^(o)t^(o) 
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C. Nonlinear terms in PRCC 

The calculations with the LPRCC approximation in- 
volves few many-body diagrams, and it is computation- 
ally less complex and hence faster. In our calculations, 
the LPRCC equations are solved first and we use the so- 
lutions as the initial guess to solve the PRCC equations. 
To describe the PRCC equations in detail, we examine 
each of the nonlinear terms. To begin with consider the 
second term on the left hand side of Eq. (21 ) and ([22|), 



second order in T, in CCSD approximation it expands to 
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All the terms contribute to both T^^'' and Similarly, 



the third term on the left hand side of Eq. (21 ) and ( 22 ), 



third order in cluster amplitude, expands to expands to 



HnT^"^T^°'^T^^^\ (24) 



In this equation, out of the four terms, only the first one 

contributes to T^^''. But, all the terms contribute to the 

T2"'^''. At the fourth order there is only one term and it 
.(1) 
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FIG. 2. Diagrams of T^^^ arising from HtiT^"'T\ 



(O)rf(l) 



1. Tj^' diagrams 

In this section we describe the single excitation dia- 
grams arising from the non linear terms. The many-body 
diagrams or the Goldstone diagrams are drawn and eval- 
uated as described in ref. [26]. The evaluation of the 
diagrams consists of calculating the radial and angular 
integra ls. C onsider the first term on the right hand side 



contributes to only T2''' ■ The terms on the right hand of Eq. (|23|), HjqTl^''T\^\ it is equivalent to ten diagrams 
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and these are shown in Fig. [2] Algebraically, we can 
write it as 



bcqa 
bpqr 

wheie gijki — (u |l/?'i2+5^(''i2)|^0 the matrix element 
of the electron-electron interactions and gijki = Qijki — 
Qijik is the antisymmtrized matrix element. We have 
used (• • • )P to represent the matrix element • • • |$o)- 
The diagrams in Fig. [2|^i-j), arising from the one-body 
part of i?Ar, evaluate to zero when orbitals are cal- 
culated with Dirac-Hartree-Fock potential. The next 

term, H^Tf^'T'"2 \ generates eight diagrams and these 
are shown in Fig. [3j It is to be noted that contractions 
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FIG. 3. Diagrams arising from the contraction _H'nTi''''T2^'' 



with only the gabpq type of two-body interaction are non- 
zero. The algebraic expression of the diagrams is 



hcqr 

^ '■b^ac ) ■ 



At the second order H^T^'^'t'^^^ is the last term. Like 
the previous term, after contraction it generates eight 
diagrams and these are shown in Fig. |4] The topological 
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(b) 



(d) 



At the third order, as mentioned earlier, only 



iJNTi^°^ri°^T^^^ contributes to the T^^^ diagrams. This 
term generate six Goldstone diagrams and these are 
shown in Fig. [5] The algebraic expression of the dia- 
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Diagrams arising from the contraction 



grams is 



E 9bcqr {tatl^ 
hcqr 



''b''a' c 



titter:). 



In total, the nonlinear terms in the T^^'' cluster equation 
generate thirty Goldstone diagrams. Considering that 



ri°) and T^^ are the dominant cluster operators in the 
unperturbed RCC and PRCC, respectively, we can ex- 
pect the mag nitude of Hj^T^°^t\^^ to be largest. 
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diagrams 



In this section we discuss the Goldstone diagrams of 
arising from the non-linear terms on the left hand 
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side of Eq. (22|. Consider the second order term, after 



expansion there arc four terms as given in Eq. ( 23 1 and 



all have nonzero contribution to T^^\ The first term. 
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FIG. 4. Diagrams arising from the contraction HjqT2'^^ T^^^ 



FIG. 6. Diagrams arising from the contraction H-f^T^"^ T^^^^ 



structure of the diagrams are very similar to those of Fig. 
[3j The algebraic expression of the diagrams is 



\nnJ^2 -"-1 /g — / , ybcqr y'^ba'p + '-fce^a + ''ab^c 
bcqr 



H^T^^'^ 't'^^\ has six diagrams and these are shown in 
Fig. [6j The equivalent algebraic expression is 



cd 



)a? represent the matrix ele- 



where, we have used 

ment ($^^| • • • |$o)- The next term, HnT^°H^^\ has six- 
teen diagrams and these are shown in Fig. [Tj However, 
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FIG. 7. Diagrams arising from the contraction Hi^T^ T 
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FIG. 8. Diagrams arising from the contraction Hj^iT.^^^ T^^^ 
^""^A/ X^M^ v^f^l^ 



the last two diagrams are zero when Dirac-Hartree-Fock- 
Brcit orbitals are used, like we do in the present work, 
in the PRCC calculations. The equivalent algebraic ex- 
pression is 



(^N^l ^2 la, 
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of T 
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^''tVca ^cVca ^ ''c^ad ^ ''a'cd I ' 

is the antisymmetrised amplitude 
Interchanging the order of excitations of the 

(oU(i) 



cluster operators, we get the next term i/N72 T^i ■ Like 
in the previous term there are sixteen diagrams and these 
are shown in Fig. [8] and equivalent algebraic expression 
is 
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,pr s 
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tabrn + T^gcdrb (CrJ 
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FIG. 9. Diagrams arising from the contraction Hj^iT^"^ 



expression for the diagrams is 
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I ^rs pq I ,pq rs\ 
+ '^abT^cd + ''cd^abl ■ 



Collecting all the diagrams, at the second order, there are 

where, tg^ = t^c " ^a? is the antisymmetrised amplitude Goldstone diagrams in the T^'^ equation after contrac- 

If., I 1 I tion of the cluster operators with iJ^r. 



of T^°\ The last second order term is HnT^°''T^^^^ and 
we can expect a large number of diagrams as both of the 
cluster operators are double excitation. There are sixteen 
diagrams and these are shown in Fig. |9] The algebraic 



tion of the cluster operators with H 

At the third order, all the terms in Eq. (24) have 



non-zero contributions to tI^-*. There are six Goldstone 



diagrams from the first term H-t^T^'^''Tl^'''T^j^^ and these 
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FIG. 10. Diagrams arising from the contraction 



are shown in Fig. |10[ The equivalent algebraic expression 
of the diagrams is 
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FIG. 12. 



Diagrams arising from the contraction 



The overall contribution from these diagrams is expected 
to be small as these are quadratic in T^'^ . The next term, 

H^T[°hl°^¥^\ has ten Goldstone diagrams and these 
are shown in Fig. |11[ The equivalent algebraic expression 



/ TT Tn(0)rTi(0)rp(l)\pg 
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-i2< - 
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'^ab^d + ^ 
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Among the third order terms in the T^2^ equation this 



is) 



will be the leading term as it depends on Tj^-* and T^^ 
the dominant among the unperturbed and perturbed 
cluster operators, respectively. There are two Goldstone 
diagrams from the fourth order term and these are shown 
in Fig. [13] and the algebraic expression is 
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FIG. 11. 



of the diagrams is 



Diagrams arising from the contraction 
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FIG. 13. 



~^'^di^a'''cb + ^c'''ab) ~ (^c^a ~ ^c^a)''"d6 

-itx-rxyj]. 

Contributions from these diagrams will be lower than the 
previous set as these depend on T2^'', which is smaller 



in magnitude than T^i \ and quadratic in T^"' like the 

previous set. The last third order term, H^T^"^ T^"^t[^\ 
has eighteen diagrams and these are shown in Fig. |12[ 
The algebraic equivalent of these diagrams is 
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Diagrams arising from the contraction 
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Among all the diagrams considered so far these two di- 
agrams will have the lowest contributions as these are 
third order in T^^\ However, for completeness we in- 
clude these in the calculations. 
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3. Dr(0' and Dr(«'r(0' dtag rams 

Another group of T^^^ diagrams arise from the contrac- 
tion of D and T^°\ which are present on the right hand 

(1) 



side of Eq. (|2lj) and ([22|). In this group, for the 
equation there are five Goldstone diagrams and these are 
shown in Fig. 14 However, among these only the last 
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FIG. 14. Singles diagrams arising from the contraction 
i^i^^fo) and i/intr('"r(o' 

one is nonlinear in T^^"^ . The algebraic expression of the 
diagrams is 



q c 
bq 

where, Tij ^ (i|rb ) is the electronic part of the single par- 
ticle matrix element. For T2^'', there are four diagrams 
and these are shown in Fig. 15 and last two are nonlinear 
in T^^"^ . The algebraic expression of the diagrams is 







(a) (b) (c) (d) 

FIG. 15. Doubles diagrams arising from the contraction 
C^rfO' and ffintTC'rC) 
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'■cry '-a'-cb ''c'^abJ' 



This completes the diagrammatic and algebraic analysis 
of the nonlinear terms in the T^-*^) equations. To obtain 
the linear algebraic equations of the cluster amplitudes, 
each of the diagrams or terms in the algebraic expres- 
sion requires further simplification to radial and angular 



components. The angular part is evaluated diagrammat- 
ically but the diagrams are different from the Goldstone 
diagrams. 



IV. DIPOLE POLARIZABILITY 

From the second order time-independent perturbation 
theory, the ground state dipole polarizability of a closed- 
shell atom is 



(^o|D|^j)(^,|D|7Ao) 

Eq — Ej 



(27) 



where |^/) are the intermediate atomic states and Ei is 
the energy of the atomic state. As D is an odd parity 
operator, l^*/) must be opposite in parity to |5'o)- In the 
PRCC theory we can write 



(^o|D|^o 
(«'o|*o) 



(28) 



From the definition of 15*0) in Eq. ^ and based on 
the parity selection rules, only the terms linear in A are 
nonzero. That is. 



(*o|*o) 
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where, D = e De , represents the unitary trans- 
formed electric dipole operator. Retaining the the lead- 
ing order terms, we obtain 
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where TV = ($o| exp[r(o)t] exp[r(o)] |$o) is the normal- 
ization factor, which involves a non-terminating series of 

contractions between T^^^^ and T'^^\ However, in the 

present work we use J\f w ($o|Ti°^^Tf ' + T^°'''' T^°'' \%) ■ 
From the above expression of a, an evident advantage 
of calculation using PRCC theory is the absence of sum- 
mation over I'^i). The summation is subsumed in the 
evaluation of the T*^^) in a natural way. This is one of 
the key calculational advantage of using PRCC theory. 

TABLE I. Comparison between GTO and GRASP92 



Atom 


GTO 


GRASP92 


Ar 


-528.6837 


-528.6837 


Kr 


-2789.8605 


-2788.8605 


Xe 


-7446.8976 


-7446.8976 


Rn 


-23602.0202 


-23602.0232 
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V. CALCULATIONAL DETAILS 

A. Basis set and nuclear density 

The first step of our calculations, which is also true of 
any atomic and molecular calculations, is to generate an 
orbital basis set. For the present work, we use the Dirac- 
Hartree-Fock Hamiltonian and even-tempered Gaussian 
type orbitals (GTOs) |17j- The radial part of the spin- 
orbitals, the large component in particular, are linear 
combinations of the Gaussian type functions 



(31) 



where p is the GTO index and C^^ is the normalization 
constant. The exponent depends on two parame- 
ters olq and /3, these are related as = ao/3^~^, where 
p = 0, 1 . . . m and m is the number of the gaussian type 
functions. The small components of the spin-orbitals are 
linear combination of g^p{r), which are generated from 
gf^p{r) through the kinetic balance condition [2H]- We 
calculate the GTOs on a grid [33] and optimize the val- 
ues of ao and /3 for individual atoms to match the spin- 
orbital energies and self consistent field (SCF) energy of 
GRASP92 [5D], which solves Dirac-Hartree-Fock equa- 
tions numerically. The comparison of the SCF energies 
are given in Table. [l] Except for Rn, there is excellent 
agreement between the SCF energies obtained from GTO 
and GRASP92. The symmetry wise values of the opti- 
mized ao and /3 are listed in Table. |llj To optimize the 

TABLE II. The ao and /3 parameters of the even tempered 
GTO basis used in the present calculations. 



Atom 


s 




P 




d 








P 


ao 




ao 


P 


Ar 


0.00055 


1.620 


0.00515 


2.405 


0.00570 


2.850 


Kr 


0.00015 


2.015 


0.00945 


2.975 


0.00635 


2.845 


Xe 


0.00012 


2.215 


0.00495 


2.995 


0.00745 


2.460 


Rn 


0.00010 


2.280 


0.00671 


2.980 


0.00715 


2.720 



basis set size, we examine the convergence of a using the 
LPRCC theory. We start with a basis set of 50 GTOs 
and increase the basis set size in steps through a series 
of calculations. As an example the results for the case of 
Kr is listed in Table. |III| The value of a changes by only 
7x 10~^ when the number of basis states is increased from 
117 to 131. So, we can use the former for our calculations 
without compromising the desired accuracy. 

In the present work we have considered finite size Fermi 
density distribution of the nucleus 



PO 



1 + e^^^^')/'^ ' 



(32) 



where, a = t41n(3). The parameter c is the half charge 
radius so that Pnuc(c) — Pq/2 and t is the skin thickness. 
The PRCC equations are solved iteratively using Jacobi 
method, we have chosen this method as it is easily paral- 
lelizable. The method, however, is slow to converge. So, 



TABLE III. Convergence pattern of a (Kr) as a function of 
the basis set size. 

No. of orbitals Basis size a 



79 


(15s,9p,9d,7/,7g) 16.8759 


97 


(17s, lip, lid, 9/, 9g) 16.7507 


117 


(21s,13p, 13d,ll/, llg) 16.7403 


131 


(25s,15p, 14d,13/,ll5) 16.7396 


139 


(25s, 16p, 15d, 13/, 135) 16.7394 


155 


(29s, 17p, 16d, 15/, 155) 16.7394 



we use direct inversion in the iterated subspace (DIIS) 
[3T] to accelerate the convergence. 



B. Breit interaction 

There are two different but equivalent approaches, re- 
ported in previous works, to calculate the matrix ele- 
ments of (7^(ri2). The first approach [32] is to couple 
the angular part of the orbitals with a to give a linear 
combination of vector spherical harmonics. This is then 
combined with the angular part of l/ri2 for integration. 
In the second approach p53], g^{ri2) is expanded as a lin- 
ear combination of irreducible tensor operators. In the 
present work we use the later and employ the expressions 
given in ref. [34] to incorporate g^{r 12) in the GTO and 
RCC calculations. For GTO calculation [351135] provides 
a very good description to include g^{r 12 in finite ba- 
sis set calculations. To assess the relative importance 
of Breit interaction, we calculate the first order energy 
correction 



IDF 



(*o| 



i<j 



(33) 



where, |$o) is the ground state reference function gen- 
erated from the Dirac-Hartree-Fock spin-orbitals and 
= "^i^j g^i^ij) represents the many-particle form 
of the Breit interaction. The (i/^)DF of the rare gas 
atoms Ar, Kr, Xe and Rn are listed in Table. IV For 



TABLE IV. SCF Energies for noble gas atoms 



Atom 


-C'SCF 






Ref. [37] 


Ar 


-528.6837 


-528.5511 


0.1326 


0.1324 


Kr 


-2788.8605 


-2787.4310 


1.4295 


1.4268 


Xe 


-7446.8976 


-7441.1248 


5.7728 


5.7753 


Rn 


-23602.0202 


-23572.8480 


29.1722 


29.3968 



each atom we calculated the SCF energy with H^'~^ and 
ffDCB^ these are E^§p = ($o|i?°^|*o) and E^§^ = 
($o|i?°™|$o)- Here, H^'^ = 7f°CB _ ^^le atomic 
Hamiltonian without the Breit interaction. From the ta- 
ble, it is evident that our results are in very good agree- 
ment with the previous results [37]. The largest devi- 
ation is observed in Rn, our result of {H^)bf is 0.8% 
lower than the previous result. However, as the Breit in- 
teraction contribution to £^gcj? is a mere 0.12% in Rn, 
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in absolute terms, the deviation is ~ 0.001%. Our re- 
sults are also in good agreement with the results of an- 
other previous study [35] ■ In the PRCC calculations, as 
described earlier, we treat at par with the residual 
Coulomb interaction. However, to examine the relative 
importance of Breit interaction, we calculate a with and 
without H^. 



VI. RESULTS AND DISCUSSIONS 

The expression of a in PRCC theory, as mentioned 
earlier, is a non-terminating series. However, the terms 
of order higher than quadratic in T have negligible con- 
tributions. For this reason, in the present work, we con- 
sider upto second order in cluster operator. As men- 
tioned earlier, we have added Breit interaction in the to- 
tal atomic Hamiltonian, one immediate outcome is, the 
number of two-electron integrals is larger and storing 
these, for faster computations, require larger memory. 
At the first order MBPT, which we use as the initial 
guess, there is an important change with the inclusion of 
H^. With only the Coulomb interaction, at the first or- 
der MBPT, the wave operator follows the parity selection 
rule and only selected multipoles of the Coulomb interac- 
tion contributes. However, with , which has opposite 
parity selection rule compared to Coulomb interaction, 
all multipoles of the two-electron interaction which sat- 
isfy the triangular conditions are allowed. In table |V) 
we list the values of a calculated using the LPRCC the- 
ory. For comparison we have also included the results 
from previous theoretical studies and experimental data. 
There are no discernible trends in the previous theoreti- 
cal results and present work. For Kr and Xe, the results 
from the many-body perturbation theory (MBPT) [T5] 
is higher than the experimental data, but with RCCSD 
triples (RCCSDT) approximations [20], Ar and Kr have 
higher values. For Ar our result is 1% higher than the ex- 
perimental data and this is consistent with the RCCSDT 
result reported in a previous work. It must, however, be 
mentioned that the previous work is based on third-order 
Douglas-KroU [21] method. Our result for Kr is in excel- 
lent agreement with the experimental data. This could 
be a coincidence arising from well chosen basis set pa- 
rameters and inherent property of PRCC to incorporate 
correlation effects more completely within a basis set. In 



TABLE V. The static dipole polarizability, a (atomic units), 
from linearized PRCC and comparison with previous results 



Method 


Ar 


Kr 


Xe 


Rn 


RCCSDTtSg 


11.22 


16.80 


27.06 


33.18 


CCSDT [H] 


11.084 


16.839 


27.293 


34.43 


MBPT[I8] 


11.062 


17.214 


28.223 




This work 


11.213 


16.736 


26.432 


35.391 


Expt.fSR] 


11.091 


16.740 


27.340 




Expt. 40^ 


11.081(5) 


16.766(8) 







mental data and 2.4 % lower than the RCCSDT result. 
The later, difference from the the RCCSDT result, can 
be partly attributed to the triple excitations. There is 
no experimental data of a for Rn, the highest Z atom 
among the noble gases. In ref. [20J, the a of Rn is com- 
puted using RCCSDT and their result is 6.2% lower than 
our result. 

To estimate the importance of Breit interaction, we 
exclude in the PRCC calculations and the value after 
excluding are 11.202, 16.728, 26.404, 35.266 for Ar, 
Kr, Xe and Rn respectively. These represent a decrease 
of 0.010, 0.012, 0.021 and 0.133 from the resuhs with 
the exclusion of H^. Except for Rn, the change in a is 
below 0.1%. This implies that to obtain accurate results 
for Rn, it is desirable to include Breit interaction in the 
calculations. 

TABLE VI. Contribution to a from different terms of the 
dressed dipole operator in the linearized PRCC theory 



Contributions from Ar 



Kr 



Xe 



Rn 



t[^^^'D + h.c. 


12.191 


18.613 


30.855 


41.641 


T,(i)tDTW + h.c. 


-0.545 


-0.888 


-1.677 


-2.328 


T2(i)tDTi°' + h.c. 


0.510 


0.748 


1.352 


1.862 


T^(l)tDy(0) + h.c. 


-0.057 


-0.118 


-0.357 


-0.301 


T2(i)tDT(°' + h.c. 


0.022 


0.038 


0.092 


0.073 


Normalization 


1.081 


1.099 


1.145 


1.157 


Total 


11.213 


16.736 


26.432 


35.391 



To examine the results in more detail, the contribu- 
tions from the terms i n th e expression of a given in Eq 
( 30 1 are listed in Table 



VI 



It is evident that T^^''^D and 
it's hermitian conjugate have the leading order contribu- 
tions. This is to be expected as these terms include the 
Dirac-Hartree-Fock-Breit contribution and RPA effects, 
which have the dominant contributions. In all the cases, 
the result from tJ^''^D is larger than the total value and 
shows dependence on Z: the results of Ar, Kr, Xe and Rn 
from this term are 8.7%, 11.2%, 16.7% and 17.7% higher 
than the total values, respectively. The next to leading 
order terms are Ti^^^^DT2'^'' and its hermitian conjugate. 
Contributions from these terms are, approximately, a fac- 
tor of twenty smaller than the leading order terms and 
opposite in phase. On closer inspection, it is natural 
that Ti(i)tDr2^°) and it's hermitian conjugate are the 
next to leading order terms. Among the second order 
terms, these are the ones which have Ti^^^ and Tj"'', the 
dominant cluster amplitudes in the perturbed and unper- 
turbed relativistic coupled-cluster theories. The results 
from Ti'^^^'I'DTj^'' have large cancellations with the term 

T2^^^^DT2°'', which is almost the same in magnitude but 
opposite in sign. Interestingly, a similar pattern is also 

2 terms. Namely, the re- 
are negative and opposite in sign 



observed with the T^i^I'Dt/"' 
suits from Ti(i)tDT2 



the case of Xe our result is 3.4% lower than the experi- 



(0) 
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TABLE VII. Contribution to a from different terms of the 
dressed dipole operator in the non-hnear PRCC theory 



Contributions from 



Ar 



Kr 



Xe 



Ti^^'I'D + h.c. 


12.950 


18.622 


33.108 


T^CDtoyOT + h.c. 


-0.579 


-0.899 


-1.7964 


T2(i)tDrf '} + h.c. 


0.488 


0.769 


1.278 


T^{i)tDr{0) ^ h.c. 


-0.061 


-0.096 


-0.392 


T2{i)tDrf ' + h.c. 


0.022 


0.035 


0.095 


Normahzation 


1.081 


1.099 


1.145 


Total 


11.859 


16.771 


28.203 



The results from the full PRCC, including the terms 



nonlinear in cluster amplitudes are given in table VII 
From the table, it is clear that the nonlinear terms tend 
to increase the deviations from the experimental data. 
A similar trend was reported in our previous work on 
Ne [H]. For Ar, the non-linear PRCC theory result is 
5.4% larger than the result from linearized PRCC and 
it is 6.5% larger than the experimental result. Similarly, 
for Xe the nonlinear PRCC result is 6.3% larger than 
the linearized PRCC result. On the other hand for Kr, 
the non-linear PRCC results are marginally larger than 
the linearized PRCC results. The larger values of a in 
the non-linear PRCC can almost entirely be attributed 
to higher value of T^^^'^D and it's hermitian conjugate. 
It means that the non-linear terms tend to increase the 
RPA effects. This is an example where inclusion of higher 
order terms enhance the uncertainty. It is possible that 
triple excitations, higher order excitation not considered 
in the present work, may balance the deviations and bring 
the results closer to the experimental data. 



addition, as the matrix elements in the expression of a 
has a quadratic dependence on radial distance, orbitals 
with larger radial extent have higher contributions. The 
next largest values arise from the npi/2 valence orbitals. 
An interesting pattern is to be noticed in the results, 
with higher Z the ratio of the contribution from 
to increases. For Ar, Kr and Xe the ratios are 2.1, 

2.3 and 2.7, respectively. However, the ratio for Rn is 
much larger, it is 5.3. The reason for the trend in the 
ratios is the contraction of the npi/2 core orbitals due to 
relativistic corrections. Hence, the npi/2 valence orbitals 
of higher Z atoms show larger contraction and accounts 
for the higher ratio. The third largest contributions in 
Ar and Kr arise from the 3si/2 and 4si/2 orbitals, respec- 
tively. This is expected as these are the orbitals which 
are energetically just below the np orbitals and spatially 
as well. Extending the same pattern, for Xe and Rn, the 
third largest contributions must be from the 5s 1/2 and 
6si/2 orbitals, respectively, but this is not case. These 
orbitals are contracted because of relativistic corrections 
and the diffused ndc, /2 orbitals have the third largest val- 
ues. From the trends in the results of the RPA effects, it 
is obvious that the relativistic corrections are important 
for Xe and Rn. 



TABLE IX. Core orbitals contribution from Ti'^'^^^DT. 
a of Argon and Krypton 



(0) 



to 



Ar 


Kr 


-0.124 (3p3/2,3pi/2) 

-0.118 (3p3/2,3p3/2) 

-0.027 (3pi/2,3pi/2) 
-0.006 (3p3/2,3si/2) 


-0.205 (4p3/2,4pi/2) 

-0.193 (4p3/2,4p3/2) 

-0.038 (4pi/2,4pi/2) 

-0.008 (4p3/2,3d5/2) 



TABLE VIII. Core orbital contribution from T^^'^D to a 



Ar 


Kr 


Xc 


Rn 


8.152 (3p3/2) 
3.914 (3pi/2) 
0.100 (3si/2) 
0.012 (2p3/2) 


12.872 (4p3/2) 
5.572 (4pi/2) 
0.058 (4si/2) 
0.056 (3d5/2) 


22.292 (5p3/2) 
8.120 (5pi/2) 
0.222 (4d5/2) 
0.140 (4d3/2) 


34.524 (6P3/2) 
6.502 (6P1/2) 
0.382 (5d5/2) 
0.214 (5d3/2) 



For a more detailed analysis of the contributions from 
the RPA effects, we consider contributions from each of 
the core orbitals in Ti^-'^D. In terms of orbital indices 



the expression is 
TWtD ^ 



H.c 



ap ' a 



(34) 



ap 



where, 



dap — (a|d|p) with d as the single particle elec- 



tric dipole operator. The four leading core orbitals (a) 
for each of the atoms are listed in Table. IVIIII In all 
the cases, the result from the outermost np2/2 valence 
orbitals are the largest. This is not surprising as these 
are the orbitals which are spatially most extended. In 



Next, we examine the pair-correlation effects, which 
manifest through the next to leading order terms, 

Ti(^)^^DT2'^'' and it's hermitian conjugate. In terms of 
orbital indices 



T^(i)tDTi") 



H.c. 



abpq 



4- fP'^ 



irPd,t-r^dp,)]. (35) 



The results of the four leading terms, listed in terms 
of the pairs of the core orbitals (ab), for Ar and Kr 
are given in Table. 



IX 



From the table we can iden- 
tify {np3f2,npi/2) as the most dominant pairing of the 
core-orbitals among the double excitations. Considering 
that the pairing is between different orbitals, the number 
of cluster amplitudes is large and this explains the large 
contribution. The second and third dominant contribu- 
tions, from the (^^3/2,71^3/2) and (npi/2, npi/2) pairs, 
are also on account of number of cluster amplitudes. 
Since ^^3/2 and npi/2 each accommodate four and two 
electrons each, respectively, the former has a larger num- 
ber of cluster amplitudes. There is a small but important 
change in the results of Xe and Rn listed in Table. |X) The 
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most dominant pair for these atoms is (np3/2, 'T-P3/2) and 
('^P3/2i is the second. This is in contrast to the 
sequence observed in Ar and Kr. The reason is, although 
the later pair has more cluster amplitudes, the npi/2 is 
contracted due to relativistic corrections. So, the contri- 
butions to a from Tj"'' involving npi/2 is smaller. The 
difference between the results from the two pairs is even 
more dramatic in Rn. There are other changes in the 
case of Rn. The (6^3/2, 5c?5/2) pair, involving the dif- 
fused 5(i5/2, is now the third largest contribution. And 
the (6pi/2, 6^1/2)7 which has the contracted orbital, 
is the fourth largest contribution. This difference in the 
sequence of leading contributions for Rn arises from the 
larger relativistic corrections. 



TABLE X. Core orbitals contribution from Ti'^'tDy^CO)-}. 
a of Xenon and Radon 



Xe 


Rn 


-0.361 (5p3/2,5p3/2) 

-0.359 (5p3/2, 5pi/2) 
-0.054 (5pi/2,5pi/2) 
-0.035 (5p3/2,4d5/2) 


-0.591 (6p3/2,6p3/2) 

-0.387 (6p3/2,6pi/2) 

-0.071 (6p3/2,5d5/2) 

-0.036 (6pi/2,6pi/2) 



To estimate the uncertainty in our calculations, we 
have identified few important sources of uncertainty. The 
first one is the truncation of orbital basis sets. Although 
we start with 9 symmetry for all the calculations, we in- 
crease the number of symmetries upto 13 in steps. The 
basis set chosen for the results given are after the value 
of a converges to 10""*. So, the uncertainty from the 
basis set truncation is negligible. The second source of 
uncertainty is the truncation of CC theory at the single 
and double excitation for both the unperturbed and per- 
turbed RCC theory. Based on earlier studies, the contri- 
butions from the triples and quadruple excitations could 
be at the most «3.3%. This is also consistent with the de- 
viations from the experimental data. Finally, the trunca- 
tion of e"^' ' De"^'"' -|-e'^' ' De"^' ' is another source of un- 
certainty. From our earlier studies with iterative method 
[42], to incorporate higher order terms in the properties 
calculations with CC theory, the contributions from the 



third or higher order in is negligibly small. Quantum 
electrodynamical (QED) corrections in this set of cal- 
culations is another source of uncertainty. However, it is 
expected to be smaller then the correction from the Breit 
interaction. As the largest Breit correction, in the case 
of Rn, is 0.1%, we can assume the corrections from QED 
effects to be at the most 0.1%. So, adding this the maxi- 
mum uncertainty in our calculations is 3.4%. However, it 
must be emphasized that, for Ar and Kr, the uncertainty 
is much smaller than this bound. 

VII. CONCLUSION 

The PRCC theory is a general extension of the RCC 
method to incorporate an additional perturbation. The 
present work demonstrates that it is suitable for prop- 
erties calculations for closed-shell atoms. Although, in 
the present work we have used PRCC to calculate elec- 
tric dipole polarizability, the method can be extended to 
calculate other atomic properties as well. 

From the present study, through the detailed analysis 
and identification of the dominant contributions, the rel- 
ativistic correction to a of noble gas atoms arising from 
the contraction of the outermost pi/2 is significant. The 
notable impact of this is the higher fractional contribu- 
tion from np^/2 in the terms which subsume RPA effects, 

T^^-'^D and it's hermitian conjugate, as we go from Ar 
to Rn. For Rn, the effect of relativistic corrections is also 
identifiable without ambiguity in the pair-correlation ef- 
fects, the (6pi/2, 6P1/2) pair is below the (6^3/2,5^5/2) 
pair for Ti'^^^^T>T!f\ 

We have also examined the importance of Breit inter- 
action in the calculation of a. The largest change of 0.1% 
is associated with Rn, the heaviest noble gas atom. So, 
when the required uncertainty of the calculations is be- 
low 1%, the inclusion of Breit interaction is desirable for 
higher Z closed-shell atoms like Rn. 
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